Horizon dynamics of distorted rotating black holes 
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We present numerical simulations of a rotating black hole distorted by a pulse of ingoing gravita- 
tional radiation. For strong pulses, we find up to five concentric marginally outer trapped surfaces. 
These trapped surfaces appear and disappear in pairs, so that the total number of such surfaces at 
any given time is odd. The world tubes traced out by the marginally outer trapped surfaces are 
found to be spacelike during the highly dynamical regime, approaching a null hypersurface at early 
and late times. We analyze the structure of these marginally trapped tubes in the context of the 
dynamical horizon formalism, computing the expansion of outgoing and incoming null geodesies, as 
well as evaluating the dynamical horizon flux law and the angular momentum flux law. Finally, we 
compute the event horizon. The event horizon is well-behaved and approaches the apparent horizon 
before and after the highly dynamical regime. No new generators enter the event horizon during 
the simulation. 



I. INTRODUCTION 

In the efforts by the numerical relativity community 
leading up to the successful simulation of the inspiral 
and merger of two black holes, analyses of single black 
holes distorted by gravitational radiation have offered a 
convenient and simpler setting to understand the non- 
linear dynamics during the late stages of binary black 
hole coalescence. For this purpose, initial data for a 
Schwarzschild black hole plus a Brill wave was presented 
in which was both time symmetric and axisymmet- 
ric. In highly distorted cases, the apparent horizon could 
develop very long, spindlelike geometries. If the event 
horizon can show similar behavior, this would raise in- 
triguing questions related to the hoop conjecture 0. The 
work of [l[ was extended to distorted rotating black holes 
in Q , where the apparent horizon served as a useful tool 
to examine the quasinormal oscillations of the black hole 
geometry as it relaxed in an evolution. Further studies 
have extracted the gravitational waves emitted by the 
black hole [3| , and compared the apparent and event hori- 
zons p|. 

We continue this line of investigation here, while incor- 
porating various modern notions of quasilocal horizons 
that have emerged in recent years. Our emphasis is on 
horizon properties during the highly dynamical regime, 
and no symmetries are present in our initial data and evo- 
lutions. The utility of quasilocal horizons can be immedi- 
ately appreciated when one wants to perform a numerical 
evolution of a black hole spacetime. One must be able 
to determine the surface of the black hole at each time, 
in order to track the black hole's motion and compute 
its properties, such as its mass and angular momentum. 
However, the event horizon, which is the traditional no- 
tion of a black hole surface, can only be found after the 
entire future history of the spacetime is known. 

Quasilocal horizons can be computed locally in time, 
and so are used instead to locate a black hole during the 
evolution. Of particular interest is a marginally outer 
trapped surface (MOTS), which is a spatial surface on 



which the expansion of its outgoing null normal van- 
ishes @. The use of MOTSs is motivated by several 
results. When certain positive energy conditions are sat- 
isfied, an MOTS is either inside of or coincides with an 
event horizon 0, 0] ■ The presence of an MOTS also im- 
plies the existence of a spacetime singularity §. In an 
evolution, the MOTSs located at successive times foliate 
a world tube, called a marginally trapped tube (MTT). 
MTTs have been studied in the context of trapping hori- 
zons [t| fToft isolated horizons [HHIIJ, and dynamical 
horizons |14| - |l6j . 

Both the event horizon and an MTT react to infalling 
matter and radiation, although their behaviors can be 
quite different in highly dynamical situations. Being a 
null surface, the evolution of the event horizon is gov- 
erned by the null Raychaudhuri equation [T3|, so that 
even though its area never decreases, in the presence 
of infalling matter and radiation the rate of growth of 
its area decreases and can even become very close to 
zero [l8|. Since an MTT is determined by quasilocal 
properties of the spacetime, its reaction to infalling mat- 
ter and radiation is often much more intuitive. A MTT 
is usually spacelike (e.g. a dynamical horizon) in such 
situations, although further scrutiny has revealed that 
MTTs can exhibit various intriguing properties of their 
own. For example, an MTT may become timelike and 
decrease in area [19|, or even have sections that are par- 
tially spacelike and partially timelike [20j | . In a numer- 
ical simulation, such behavior is often indicated by the 
appearance of a pair of new MTTs at a given time, ac- 
companied by a discontinuous jump in the world tube of 
the apparent horizon, or outermost MOTS. 

In this paper, we investigate the behavior of MTTs 
and the event horizon in the context of a rotating black 
hole distorted by an ingoing pulse of gravitational waves. 
First, we construct a series of initial data sets in which the 
amplitude of the gravitational waves varies from small to 
large, which are then evolved. We focus on the evolution 
with the largest distortion of the black hole, in which the 
mass of the final black hole is more than double its initial 
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value. During the evolution, the world tube of the appar- 
ent horizon jumps discontinuously when the gravitational 
waves hit the black hole, and as many as five MTTs are 
found at the same time. Some of these MTTs decrease 
in area with time, although we find that all the MTTs 
during the dynamical stages of our evolution are space- 
like and dynamical horizons. Moreover, all these MTTs 
join together as a single dynamical horizon. Their prop- 
erties are further analyzed using the dynamical horizon 
flux law [l5j . which allows one to interpret the growth 
of the black hole in terms of separate contributions. We 
also evaluate the angular momentum flux law based on 
the generalized Damour-Navier-Stokes equation [2l[ ■ Fi- 
nally, we locate the event horizon and contrast its behav- 
ior with that of the MTTs. 

The organization of this paper is as follows. Section II 
details the construction of the initial data sets and Sec. 
Ill describes the evolutions. Section IV introduces some 
definitions about MOTSs, and the methods used to lo- 
cate them. Section V discusses the MTTs foliated by 
the MOTSs, the determination of their signatures, and 
the fluxes of energy and angular momentum across them. 
The emphasis is on the case with the largest distortion 
of the initial black hole, as is the remainder of the pa- 
per. Section VI explains how we find the event horizon, 
and contrasts its properties with the MTTs. Section VII 
presents some concluding remarks. Finally, the appendix 
offers some insight on our results in light of the Vaidya 
spacetime. 



II. INITIAL DATA 

Initial data sets are constructed following the method 
of [22| . which is based on the extended conformal thin 
sandwich formalism. First, the 3+1 decomposition of 
the spacetime metric is given by [23|, [24[ 

( 4 W = g^dx^dx", (1) 
= -N 2 dt 2 + g i:j (dx l + tfdt) (dx j + I3 3 dt) , (2) 

where gij is the spatial metric of a t = constant hypcrsur- 
face Tit, N is the lapse function, and j3 l is the shift vec- 
tor. (Here and throughout this paper, Greek indices are 
spacetime indices running from to 3, while Latin indices 
are spatial indices running from 1 to 3.) Einstein's equa- 
tions (here with vanishing stress-energy tensor = 0) 
then become a set of evolution equations, 



{d t -Cp) 9lJ =-2NK ih 
(d t - C f3 )K tJ = N (Rtj - 2K ik K k j + KKi 



and a set of constraint equations 

R + K 2 ~K t] K lJ = 0, 
Vj (K ij -g lJ K) = 0. 



(3) 
(4) 



(5) 
(6) 



In the above, C is the Lie derivative, Vj is the covariant 
derivative compatible with g^, R = g tJ Rij is the trace of 
the Ricci tensor Rij of gij, and K = g^ is the trace 
of the extrinsic curvature of T, t . 

Next, a conformal decomposition of various quantities 
is introduced. The conformal metric gij and conformal 
factor ip are given by 



9ij = i> 4 9ij, 



(7) 



the time derivative of the conformal metric is denoted by 



dtgij, 



(8) 



and satisfies Uijg 1 -' = 0, while the conformal lapse is given 
by N = ip~ 6 N. Equations ©, ©, and the trace of © 
can then be written as 



VV - i^R - ^ 5 K 2 + ^-r A ..#j = 0j 

(9) 

2 -^V l K = 0, 



2N 



1 



,-,'■> 



2N 
1 ^ 5 



(10) 



V 2 (VV/) - (N^) i^R + ^K 2 + L^Av® 

= -V' 5 {dtK - p k d k K) . 

(11) 

In the above, Vj is the covariant derivative compatible 
with g^, R = g^ Rij is the trace of the Ricci tensor R^ 
oigij, L is the longitudinal operator, 



L/3 



13 = V7? J +V J /f -|g y 'V fc /3 fc , (12) 



and A l 3 is 



2N 

which is related to by 



K 



i,- w A l3 + - 9ij K. 



(13) 



(14) 



For given g i3 , u l3 , K, and d t K, Eqs. ©, (fTUj) . and pT|) 
are a coupled set of elliptic equations that can be solved 
for ip, N, and j3 l . From these solutions, the physical 
initial data g^ and are obtained from (O and (|14|). 
respectively. 

To construct initial data describing a Kerr black hole 
initially in equilibrium, together with an ingoing pulse 
of gravitational waves, we make the following choices for 
the free data, 



gij = g^ s + Ahij, 



iij = Adthij ~ - 9 ij9 Ad t hku 



K = K 
d t K = 0. 



KS 



(15) 
(16) 

(17) 
(18) 
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FIG. 1. Convergence of the elliptic solver for different am- 
plitudes A. Plotted is the square-sum of the Hamiltonian 
and momentum constraints, Eqs. (0 and as a function of 
numerical resolution, measured here by the number of radial 
basis functions in the spherical shell containing the gravita- 
tional waves. 



In the above, g^ s and K KS are the spatial metric and 
the trace of the extrinsic curvature in Kerr-Schild coordi- 
nates, with mass parameter A/ks = 1 an d spin parameter 
o-ks = 0.7Mks along the z-direction. The pulse of gravi- 
tational waves is denoted by hij, and is chosen to be an 
ingoing, even parity, m = 2, linearized quadrupole wave 
in a flat background as given by Teukolsky (25( (see [26| 
for the solution for all multipolcs). The explicit expres- 
sion for the spacetime metric of the waves in spherical 
coordinates is 



hij dx % dx 3 



- (Ri sin 2 9 cos 2(f) dr 2 
2i?2 sin 6 cos 6 cos 2(f>rdrd9 
2i?2 sin 6 sin 2(f>r sin 6drd(f> 
[R 3 (l + cos 2 9) cos 2(p - R t cos 2<f\ r 2 d 2 
[2(R 1 -2R 3 ) cos 9 sin 20] r 2 sin 9d9d<p 
[R 3 (l + cos 2 9) cos 2(j) + Rx cos 2 9 cos 2q 
x r 2 sin 2 8d 2 (h, 



where the radial functions are 

- F (2) 3F (1) 3F 



Ri = 3 

R 2 = - 



i?(3) 3j p(2) 6F (1) Q F 



(19) 

(20) 
(21) 
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FIG. 2. ADM energy £adm and Christodoulou mass Mi of 
the initial data sets, versus the gravitational wave amplitude 
A. The inset shows the Ricci scalar R along the x-axis. All 
quantities are given in units of the mass of the background 
Kerr-Schild metric. 



fla = - 



i?(4) 2F^ 9F< 2 ' 2LF« 21F 



and the shape of the waves is determined by 



F 



F = F(t + r) = F(x) = e -(,*-x<y) 2 /™\ 
'd n F(x)~ 

x—t-\-r 



dx n 



(22) 

(23) 
(24) 



We choose F to be a Gaussian of width w/Mks = 1-25, at 
initial radius xq/Mks = 15. The constant A in Eq. (fT5")l 
is the amplitude of the waves. We use the values A = 0.1, 
0.2, 0.3, 0.4, and 0.5, each resulting in a separate initial 
data set. 

Equations (j9)) , (| 10[) , and (|TT|) are solved with the pseu- 
dospectral elliptic solver described in [27]]. The domain 
decomposition used in the elliptic solver consists of three 
spherical shells with boundaries at radii r/Afxs = 1-5, 12, 
18, and 10 9 , so that the middle shell is centered on the 
initial location of the gravitational wave pulse. The inner 
boundary lies inside the apparent horizon and Dirichlct 
boundary conditions appropriate for the Kerr black hole 
are imposed. It should be noted that these boundary 
conditions are only strictly appropriate in the limit of 
small A and large xq, when the initial data corresponds 
to an ingoing pulse of linearized gravitational waves on an 
asymptotically flat background, with a Kerr black hole 
at the origin. As A is increased and xq is reduced, we ex- 
pect this property to remain qualitatively true, although 
these boundary conditions become physically less well 
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motivated. Nonetheless, we show below by explicit evo- 
lution that most of the energy in the pulse moves inward 
and increases the black hole mass. 

At the lowest resolution, the number of radial basis 
functions in each shell is (from inner to outer) N r = 9, 
18, and 9, and the number of angular basis functions in 
each shell is L = 5. At the highest resolution, the num- 
ber of radial basis functions in each shell is (from inner to 
outer) N r — 41, 66, and 41, and the number of angular 
basis functions in each shell is L = 21. Figure Q] shows 
the convergence of the elliptic solver. The expected expo- 
nential convergence is clearly visible. Curves for each A 
lie very nearly on top of each other, indicating that con- 
vergence is independent of the amplitude of the waves. 
We evolve the initial data sets computed at the highest 
resolution of the elliptic solver. 

We locate the apparent horizon (the outermost 
marginally outer trapped surface defined in Sec. lIV A| in 
each initial data set using the pscudospcctral flow method 
of Gundlach [28| (explained briefly in Sec. IIVB[) . and 
compute the black hole's initial quasilocal angular mo- 
mentum Ji and Christodoulou mass Mj (the subscript 
"£" denotes initial values). The quasilocal angular mo- 
mentum J is defined in Eq. (1481) . which we calculate 
with approximate Killing vectors|29l| (see also [30] ) . The 
Christodoulou mass M is given by 



M 



'Mb 



r- 



(25) 



where Mh = -J An/lQn is the Hawking or irreducible 
mass [3l| , with Ah being the area of the marginally outer 
trapped surface of interest. The main panel of Fig. [2] 
shows M and the Arnowitt-Deser-Misner (ADM) energy 
£adM) as a function of the amplitude A of each initial 
data set. The difference between -Eadm an d M is a mea- 
sure of the energy contained in the ingoing gravitational 
waves. For A > 0.4, this energy is comparable to or 
greater than M, so the black hole will become strongly 
distorted in the subsequent evolution. The inset of Fig. [2] 
shows the Ricci scalar R of gij along the x-axis at the ini- 
tial location of the gravitational wave pulse. The sharp 
features of R necessitate the use of the higher N r as la- 
beled in Fig. Q] 



III. EVOLUTIONS 

Each of the initial data sets are evolved with the Spec- 
tral Einstein Code (SpEC) described in [H,[33j]. This code 
solves a first-order representation {3~i| of the generalized 
harmonic system [3 5143 7| . The gauge freedom in the gen- 
eralized harmonic system is fixed via a freely specifiable 
gauge source function that satisfies 



(26) 
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FIG. 3. Constraint violations for the evolution with A = 0.5. 
Plotted is the L 2 norm of all constraints, normalized by the 
L 2 norm of the spatial gradients of all dynamical fields. 



equations for N and (3 l [34j . 

d t N - frdiN = -N (H t - pHi + NK) , (27) 
dtp - P k d k p = Ng l i [N {Hj + g kl T jkl ) - djN] , (28) 

so there is no loss of generality in specifying instead of 
N and /3 % , as is more commonly done. For our evolutions, 
.Hp is held fixed at its initial value. 

The decomposition of the computational domain con- 
sists of eight concentric spherical shells surrounding the 
black hole. The inner boundary of the domain is at 
r/Afxs = 1-55, inside the apparent horizon of the ini- 
tial black hole, while the outer boundary is at t/Mks = 
50. The outer boundary conditions 34,[38|,[39| are de- 
signed to prevent the influx of unphysical constraint vio- 



lations 
diation 




and undesired incoming gravitational ra- 
48| , while allowing the outgoing gravitational 
radiation to pass freely through the boundary. Intcrdo- 
main boundary conditions are enforced with a penalty 
method [iij [5(j ■ The evolutions were run on up to three 
different resolutions - low, medium, and high. For the 
low resolution, the number of radial basis functions in 
each shell is N r = 23, and the number of angular basis 
functions in each shell is L = 15. For the high resolution, 
N r = 33 and L = 21 in each shell. 

We will be mainly interested in the case where the 
gravitational waves have an amplitude A = 0.5. As a 
measure of the accuracy of this evolution, the constraints 
of the first-order generalized harmonic system are plotted 
in Fig. [3] Plotted is the L 2 norm of all constraint fields, 
normalized by the L 2 norm of the spatial gradients of the 



where = <? r M „A is the trace of the Christoffcl sym- dynamical fields (sec Eq. (71) of [34[). The L norms 
bol. In 3+1 form, the above expression gives evolution are taken over the entire computational volume. The 
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constraints increase at first, as the black hole is distorted 
by the gravitational waves. As the black hole settles down 
to equilibrium, the constraints decay and level off. The 
results presented in the following sections use data from 
the high resolution runs only. 



B. MOTS Finders 

We use two different algorithms to locate MOTSs in 
Ti t . Both algorithms expand an MOTS "height function" 
in spherical harmonics 



IV. MARGINALLY TRAPPED SURFACES 
A. Basic Definitions and Concepts 

Let S be a closed, orientablc spacelike 2-surface in S t . 
There are two linearly independent and future-directed 
outgoing and ingoing null vectors Z M and normal to 
S. We write these vectors in terms of the future-directed 
timclike unit normal to E t and the outward-directed 
spacelike unit normal s M to S as 



l» = —= (n** + sn and k^ = —= (n* - 
y/2 \/2 



normalized so that g^l^k 11 
metric on iS is 



(29) 

1. Then the induced 



(30) 
(31) 



The extrinsic curvatures of S as embedded in the full 
four-dimensional spacetime are 

K^l = q^Vxlp and K$ = q^y X k p . (32) 

The null vectors l^ 1 and k^ are tangent to a congruence 
of outgoing and ingoing null geodesies, respectively. The 
traces of the extrinsic curvatures give the congruences' 
expansions 

= g^V^ and 6 {k) = g^V^fe, (33) 
and the shears are the trace-free parts, 

^1 = ^^x^-^,0(1) and (34) 



(35) 



The geometrical interpretation of the expansion is 
the fractional rate of change of the congruence's 
cross-sectional area [13] ■ We will mainly be inter- 
ested in 2-surfaces S on which 6 m = 0, called 
marginally outer trapped surfaces (MOTSs) following the 
terminology in [2Cj. If 6m < on S, then outgoing 
null normals will be converging towards each other, as 
one expects to happen inside a black hole. If On) > 
the situation is reversed, so the condition 6m = pro- 
vides a reasonable quasilocal prescription for identifying 
the surface of a black hole. In practice, an MOTS will 
generally lie inside the event horizon, unless the black 
hole is stationary. The outermost MOTS is called the 
apparent horizon, and is used to represent the surface of 
a black hole in numerical simulations. In the next sub- 
section, we briefly describe how we locate MOTSs. 



HviotsI 



^MOTS I 

E E 



(36) 



1=0 



Our standard algorithm is the pseudospectral fast flow 
method developed by Gundlach [28| . which we use dur- 
ing the evolution. This method utilizes the fact that the 
MOTS condition 6m = results in an elliptic equation 
for rMOTs(0,4>)- The elliptic equation is solved using a 
fixed-point iteration with the flat-space Laplacian on S 2 
on the left-hand side, which is computationally inexpen- 
sive to invert given the expansion Eq. (|36p . The fixed- 
point iteration is coupled to parameterized modifications 
which allow for tuning of the method to achieve fast, 
but still reasonably robust convergence. In Gundlach's 
nomenclature, we use the N flow method, and have found 
the parameters a = 1 and f3 = 0.5 satisfactory (see [28| 
for definitions). 

Gundlach's algorithm (as well as MOTS finders based 
on flow methods in general [5l|, HH) incorporates a sign 
assumption on the surfaces near the MOTS, namely that 
6m is positive for a surface which lies somewhat outside of 
the MOTS. This assumption is satisfied for the apparent 
horizon. However, this sign assumption is not satisfied 
for some inner MOTSs in E t that we discover below. 
Therefore, these inner MOTSs are unstable fixed-points 
for Gundlach's algorithm, so that this algorithm cannot 
locate these MOTSs. 

To find these inner MOTSs, we employ an older algo- 
rithm that is based on a minimization technique [53l455| : 
The coefficients Ai rn in Eq. (|36p are determined by min- 
imizing the functional 
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(37) 



where the surface integral is over the current trial surface 
with area element y/q. This technique is insensitive to 
the sign assumption in Gundlach's method. However, it 
is much slower, especially for large £mots- 

When multiple MOTSs arc present in St, the choice of 
an initial surface determines the final surface the MOTS 
finder converges to. Therefore, both MOTS finders re- 
quire judicious choices of these initial surfaces. We typi- 
cally track MOTSs from time step to time step, and use 
the MOTS at the previous time step as an initial guess 
for the MOTS finder at the current time. 
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V. MARGINALLY TRAPPED TUBES 
A. Basic Definitions and Concepts 

During an evolution, the MOTSs found at successive 
times foliate a world tube, or a marginally trapped tube 
(MTT). The type of MTT that is foliated by a scries 
of MOTSs depends on the physical situation. A null 
MTT is an isolated horizon juHH [H, E3 if -R^J U is 
future causal, and certain quantities are time indepen- 
dent on it. An isolated horizon describes a black hole in 
equilibrium. On the other hand, a dynamical horizon 
describes a black hole that is absorbing matter or 
gravitational radiation 0, [l5j], and is physically the 
most relevant. A dynamical horizon is a spacelike 
MTT foliated by MOTSs on which 6 {k) < 0, called 
future marginally outer trapped surfaces. For a given 
slicing of spacetimc by spatial hypersurfaces S t , the foli- 
ation of a dynamical horizon by future marginally outer 
trapped surfaces on E t is unique fl6| . Since the location 
of a MOTS is a property of S t , different spacetime slic- 
ings will in general give different MTTs. Also, a timelike 
MTT is called a timelike membrane (58j . Since causal 
curves can traverse it in both inward and outward direc- 
tions, it cannot represent the surface of a black hole. 

An additional characterization of MTTs is based on 
trapping horizons @. A future outer trapping horizon is 
an MTT foliated by MOTSs that have 6 {k) < and 
£k9(i) < for some scaling of l^ and fc M . Such an 
MOTS is called a future outer trapping surface. If the 
null energy condition holds, a future outer trapping hori- 
zon is either completely null or completely timelike. It 
was shown in [59| that if CkO(i) ^ for at least one 
point on these future outer trapping surfaces, then the 
future outer trapping horizon is spacclike, or a dynamical 
horizon, in a neighborhood of the future outer trapping 
surfaces. Otherwise the future outer trapping horizon is 
null. 

Interestingly, an MTT may not fall into either of the 
categories described above, but can have sections of 
mixed signatures as demonstrated in the head-on colli- 
sion of two black holes jlJoj . At merger, a common ap- 
parent horizon appears in S t that surrounds the MOTSs 
of the individual black holes. This common horizon then 
bifurcates into outer and inner common horizons. The 
outer common horizon grows in area and is spacelike. 
However, the inner common horizon decreases in area 
and foliates an MTT that is briefly partly spacelike and 
partly timelike, before becoming a timelike membrane 
later on. 



B. Multiple MTTs 

We now discuss the MOTSs that occur during the five 
evolutions of the distorted black hole, with amplitude 
A = 0.1, 0.2, 0.3, 0.4, or 0.5 for the ingoing gravitational 
wave pulse. The MOTSs we find are indicated in Fig. @] 




10 11 12 13 14 15 16 17 18 19 

tIM. 
i 

FIG. 4. The solid curves are the Christodoulou masses M(t) 
divided by their initial values M; for the five evolutions with 
different amplitudes A = 0.1, 0.2, 0.3, 0.4, and 0.5 for the 
ingoing pulse of gravitational waves. The horizontal dotted 
lines denote the ADM energy of each data set, Badm/ 'Mi. 



by their Christodoulou masses M. Early in each simula- 
tion, M is approximately constant, and begins to increase 
when the gravitational wave hits the black hole around 
t Kb 12 Mi. The effect is more pronounced for larger A. 
The horizontal dotted lines in Fig. @] indicate the ADM 
energy of the initial data. Although we do not explicitly 
calculate the energy carried away by gravitational waves, 
we can still see that the final Christodoulou mass is close 
to .Eadm, indicating that the energy in the gravitational 
wave pulse predominantly falls into the black hole, and 
only a small fraction of this energy propagates to null 
infinity. Even for the highest amplitude case of A = 0.5, 
the final value of M is about 99.1% of the ADM energy. 
These results are as expected. However, for both A = 0.4 
and A = 0.5, a very interesting new feature arises: multi- 
ple concentric MOTSs are present at the same coordinate 
time. 

The evolution with A = 0.5 shows the multiple MOTSs 
more distinctly, hence we will focus on it in the remain- 
der of this paper. Figure [5] presents a closer look at the 
irreducible masses Mh for this case. Locating all MOTSs 
shown in Fig. [5] requires considerable care. The start- 
ing point was the output of the MOTS finder that was 
run during the evolution, using Gundlach's fast flow al- 
gorithm [281 ] . Because of the computational expense in- 
volved, the MOTS finder was not run very frequently, 
resulting in the solid circles in Fig. [5] The MOTS at the 
previous time was used as the initial guess for the current 
time, resulting in a series of MOTSs which is as continu- 
ous as possible. The curve traced out by these points has 
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FIG. 5. Irreducible mass Mh divided by its initial value Mh> 
for the evolution with A = 0.5. The solid circles are the values 
of Mh for MOTSs found during the evolution. The completed 
curve is traced out by open circles. The vertical shaded region 
indicates when five MOTSs exist at the same time. 



FIG. 6. The solid red line denotes the apparent horizon for 
the evolution with A — 0.5. The solid blue circles denote 
an erroneous "apparent horizon," which is found when the 
apparent horizon finder is run during the evolution in larger 
time intervals. The black dashed lines denotes all five MTTs 
as shown in Fig. [5] 



sharp jumps, which was the first indication of the pres- 
ence of multiple MOTSs at these times. Then to find the 
remainder of MTT3 and MTT5, an MOTS correspond- 
ing to one of these solid circles on MTT3 or MTT5 was 
used as an initial guess and the MOTS finder was also run 
more frequently. At this stage, we had completely traced 
out MTT1, MTT3, and MTT5. Next we found MTT2 
and MTT4 to be unstable fixed points for Gundlach's 
algorithm, so it was necessary to use our older MOTS 
finder based on a minimization technique [53t455| to find 
these MTTs. As an initial guess for finding an MOTS on 
MTT2 for instance, a sphere with radius equal to the av- 
erage radii of MTT1 and MTT3 sufficed. Once an MOTS 
on MTT2 was located, it was used as an initial guess for 
the MOTS finder to locate the MOTSs on neighboring 
time slices (both later and earlier). The same procedure 
was used to locate MTT4. 

After finding all the MTTs in Fig. O a clearer pic- 
ture of their structures in relation to each other emerged. 
MTT1 corresponds to the surface of the initial black hole. 
Shortly after t = UM h a new MOTS with M H /M H ,; ~ 
1.525 appears and bifurcates into two MTTs. Mh de- 
creases along MTT2, which promptly annihilates with 
MTT1, while MTT3 persists slightly longer. A similar 
process then takes place again, and MTT5 is left over as 
the surface of the final black hole, with Mh more than 
double its initial value. The vertical shaded region indi- 
cates the time interval when five MTTs exist simultane- 
ously. Notice that Mh of the apparent horizon jumps dis- 
continuously in time from the curve of MTT1 to MTT3, 
and then to MTTS. This indicates that the apparent 



horizon itself is discontinuous across these times. 

The apparent horizon is the outermost MOTS, and 
when only one MOTS is present in a black hole evolution, 
the MOTS and apparent horizon are identical. Here this 
is not the case, and Fig. [5] shows the apparent horizon in 
relation to the various MTTs. This figure also highlights 
another potential pitfall when locating MOTSs. MOTS 
finders are typically run during the evolution fairly in- 
frequently, using the MOTS from the last MOTS com- 
putation as an initial guess (to minimize computational 
cost). If this had been done for the A = 0.5 case shown 
in Figs. [5] and El the solid blue circles would have been 
obtained. Because the previously found MOTS is used 
as an initial guess, newly appearing MOTSs are generally 
missed. For instance, the solid blue circles follow MTT1 
until it disappears, instead of jumping to MTT3. There- 
fore, the output of the "apparent horizon finder" (the 
more widely used name, but technically less precise than 
"MOTS finder"), is sometimes not the apparent horizon. 

A measure of the distortion of the black hole is pro- 
vided by the intrinsic scalar curvature R of the MOTSs. 
The extrema of R is shown in Fig. [7J along with those of 
the initial apparent horizon. It is interesting to point out 
that around t = 14.25Mi, the distortion caused by the 
gravitational waves with A = 0.5 is sufficiently strong to 
produce regions of negative R. 
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FIG. 7. Extrema of the intrinsic scalar curvature R of MOTSs 
during the evolution with A = 0.5. The horizontal dotted 
lines are the values for the apparent horizon in the initial data. 
Around t = 14.25Mj, the MOTSs have regions of negative R. 



C. Dynamical Horizons 

We determine the signatures of the multiple MTTs 
during the highly dynamical period. First we compute 
On.} and CkO(i) using the null normals in Eq. (|29[) . and 
find that both quantities arc negative. So our MTTs 
are future outer trapping horizons, which must be cither 
spacelike or null, and we can immediately rule out the 
possibility of there being sections of mixed signatures. 
Figure [5] shows the extrema of 9^ along each MTT. The 
quantity Ck0(i) is evaluated from the expression [59| 

erf® = -R/2 + - d^u" + SttT^W, (38) 
where 



= -qZk X W v l } 



(39) 



is the normal fundamental form, and d^ is the covariant 
derivative compatible with q^. Figure [9] shows the ex- 
trema of C k 9{i) < along each MTT. 

Next wc compute to determine whether the 

MTTs are spacelike or null. We evaluate this using the 
null Raychaudhuri equation (l7j . 



0-jH — ^TlT^vl^l" . 



(40) 



Figure [10] shows that during the times when there are 
multiple MTTs, C t 9 {l) ^ somewhere on each MOTS. 
Thus all of the MTTs are dynamical horizons at these 
times. 



Here we also mention the extremality parameter e of a 
MTT introduced in j6fj| . In vacuum, it is given by 



— / UnUJ^ \fqd 2 x, 
4tt J s 



is 
1 

47T 



\fqd 2 x, 



(41) 
(42) 



where the integral is over an MOTS S that foliates the 
MTT. When S is axisymmctric, this can be regarded as 
the sum of the squares of all angular momentum multi- 
poles. Because a future outer trapping horizon, which 
is either spacelike or null, has Ck0(i) < 0, it is always 
subextremal (e < 1). So a timelike membrane foliated 
by future MOTSs (with 6 {k) < 0) must have C k 9 {i) > 0, 
and is superextremal (e > 1). Therefore, it was suggested 
in [6(| that an MTT's transition from being spacelike to 
timclike can be detected when e — > 1. 

Figure [TT] shows e along each MTT, and we see that 
nowhere does e — > 1, confirming that our MTTs do not 
become timelike. The value of e shows a substantial 
decrease after the distortion has left, which is not due 
to a loss of quasilocal angular momentum J (defined in 
Eq. (|48|l). but to the large gain in irreducible mass Mh. 
It may seem that e in Fig. [14] is already rather small to 
start out with, but one must recall that e depends on 
the scaling of the null normals l 1 * and . That is, we 
can define new null normals l^ = fl 11 and — k^/f, 
rcscalcd by some function / such that the normalization 



-1 is preserved. Then e will change as 



e = e + i f [2^ I d A1 ln/ + (d A1 ln/)(^ln/)] ^dd 2 x. (43) 
4tt J s 

Nevertheless, the extremality classification of the MTTs 
is invariant. 

It is known that the irreducible mass Mh of an MOTS 
must increase along a dynamical horizon |15| , so at first 
it may seem surprising that MTT2 and MTT4, with de- 
creasing Mh during the evolution, are also dynamical 
horizons. However, all these MTTs can be viewed as 
sections of a single dynamical horizon H that weaves for- 
wards and backwards in time. Then it is clear that the 
tangent vector to H along MTT2 and MTT4 points back- 
wards in time, so that Mh is actually increasing along 
H as expected. Our simple choice of holding the gauge 
source function equal to its initial value leads to a 
spacctimc foliation that interweaves H. This could be 
avoided by an alternative choice of that results in a 
single dynamical horizon that only grows in time. 

The situation here resembles an example of a Tolman- 
Bondi spacetime considered in [r9j , where multiple spher- 
ically symmetric dust shells fall into a black hole. For 
their chosen matter distribution, multiple MTTs also 
formed (up to three at the same time), which were ei- 
ther completely spacclike, or null when the matter den- 
sity vanished between successive dust shells. In our case 
the role of the matter density is replaced by the shear a^l 
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FIG. 8. Extrema of 6 (k) on each MOTS along the MTTs 
during the evolution with A = 0.5. For the time shown, 0^) < 
0. 



FIG. 9. Extrema of C k 8 {l) on each MOTS along the MTTs 
during the evolution with A — 0.5. For the time shown, 

c k e m < 0. 



due to the gravitational waves. Since this is always non- 
vanishing somewhere on the multiple MTTs that form, 
we only have dynamical horizons. 

In [lq , it was shown that for a regular dynamical hori- 
zon (which is achronal and also a future outer trapping 
horizon), no weakly trapped surface (on which 9m < 
and 9(u\ < 0) can exist in its past domain of dependence. 
This helps to explain the difficulty in locating MOTSs 
along MTT2 and MTT4 using flow methods. For exam- 
ple, consider locating an MOTS on MTT2 at t = 14. 1M; 
shown in Fig. [5] If we use a trial surface S located be- 
tween the MOTSs on MTT1 and MTT2. it must have 
6m > because it lies in the past domain of dependence 
of %. This means that S will be moved inwards when 
using flow methods, away from MTT2. If we switch to 
having S lie between the MOTSs on MTT2 and MTT3, 
then having 9m > is desired. Unfortunately, now S lies 
in the future domain of dependence of Ji, and we are no 
longer guaranteed that S is not a weakly trapped surface. 



D. Dynamical Horizon Flux Law 

The growth of a black hole in full, nonlinear general 
relativity can be described by the dynamical horizon flux 
law of Ashtekar and Krishnan [lj, [l5| , which relates the 
increase in area or mass along a dynamical horizon to 
fluxes of matter and gravitational energy across it. Here, 
we will evaluate this flux law for the dynamical horizon 
T-L that consists of the multiple MTT sections we found 
earlier, using the form given in (59|. 

To state the dynamical horizon flux law, let us specif- 
ically consider the change in the irreducible mass Mh 



along T-L. Denote an MOTS that foliates H by S v , which 
is labeled by a foliation parameter v that is constant on 
S v . Then choose a tangent vector V M to H that is normal 
to each S v , and such that 



Cyv = 1. 



This vector can be written as 



= gin _ C ' fc M ) 



(44) 



(45) 



in terms of coefficients B and C, and null normals = 
fl^ and fc M = k^/f that are rescalcd by a function / (but 
still having l^k^ = —1) so that 



C = 2 



dv 



(46) 



The dynamical horizon flux law is then 



dv 



T V- 



S7T 



8tt 



(47) 

where a^l and <D P are given by Eqs. (|34| and (|39| but in 
terms of and and r M = Bl^ + Ck^ is the normal 
vector to H. 

The first term in Eq. ([47]) involving is the energy 
flux of matter across S v , and the second term involv- 
ing a^l is a flux of gravitational energy [l5j]. The last 
term has been interpreted differently by various authors. 
The normal fundamental form lo^ (or ii„) enters into the 
definition of the quasilocal angular momentum J of a 
black hole mentioned at the end of Sec. [Til which is given 
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FIG. 10. Extrema of -&0 W on each MOTS along the MTTs 
during the evolution with A = 0.5. Near t = 14M t , Ci6(i) / 
somewhere on S v . 



FIG. 11. Extremality parameter e along the MTTs during 
the evolution with A = 0.5. For the time shown, the MTTs 
are subextremal with e < 1, indicating that the MTTs have 
no timelike sections. 



by [l5j, 



J 



1 

'8tt 



Vqd 2 



(48) 



for any choice of rotation vector field </> M on S v . Because 
of this relation, this term has been interpreted as a flux of 
rotational energy [H, [2(j. However, it has been pointed 
out in [5^] that this is unlikely, as uj^ is related to J 
itself and not its flux. Indeed, this may be illustrated 
by considering a Kerr black hole that is distorted by an 
ingoing spherically symmetric dust shell (which carries 
no angular momentum). So even though there will be no 
flux of rotational energy, the last term in Eq. (|47p will still 
be nonzero whenever C^0, which is necessarily true on a 
dynamical horizon. This last term also closely resembles 
the extremality parameter e mentioned in Sec. IV CI 

Another interpretation of the last term in Eq. (|4T|) 
has been given by Hayward [6lJ as a flux of longitudinal 
gravitational radiation, by examining the components of 
an effective gravitational radiation energy tensor in spin- 
coefficient form. At future null infinity, the outgoing lon- 
gitudinal gravitational radiation is negligible relative to 
the outgoing transverse radiation, but near the black hole 
this is generally not so. 

To evaluate the dynamical horizon flux law, we first 
construct a tangent vector A' 1 to % that connects S v in 
Tit to S v / >v in T, t i as 



X» = ± 1 



dt 



(49) 



where Xy EIX6 the coordinates of S v , and the plus sign is 
for t' > t while the minus sign is for t' < t. The latter 



occurs along MTT2 and MTT4. The spatial components 
of the tangent vector A M diverge when two MTT sections 
meet. This may be avoided by a different choice of A M , 
but here we employ the simple one described above. For 
this reason, we also consider the corresponding foliation 
parameter v along each section of H separately. Since 



£ x v 



dv 

at' 



(50) 



and we would like this to be unity, it follows that v = ±t+ 
i>o, where vq is some constant along each MTT section. 
We choose v = t along MTT1. Along the other MTT 
sections, we choose vq so that v — on the first S v we 
find on those sections. 

Next we make A M orthogonal to S v to obtain (while 
leaving the time component unchanged, so Eq. (|44"]) is 
still satisfied with the choice of v described above). To 
achieve this, we use the unit tangent vectors to S v , 



P' 



N v 



dxl 
06 



and q>* = N„ [ 0, 



1 dxl 



O0 



sin ( 



(51) 



Here, x\{6,<j>) = c^ots + r MOTs(0, <£K(0, </>) where 
?>iots(0, <t>) is given in Eq. (|5o) and cf is the coordi- 
nate unit vector pointing from the origin c^ots °f ^ ne 
expansion along the (6, </>)-directions. Also 
are normalization factors such that p 2 = q 2 
onalizing against gives the vector 

Q*=N Q {<f-irq v j?) l (52) 

where Nq is again a normalization factor such that Q 2 = 
1. Then we obtain the desired tangent vector to H as 



N p and N q 
1. Orthog- 



V 



X»-{p v X u )p»-{Q v X v )Q> 



(53) 
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FIG. 12. Terms in the dynamical horizon flux law of Eq. (|47[) 
plotted against the foliation parameter v along each section 
of S v . Along MTT1, we choose v = t. Along the other MTT 
sections, we choose v = on the first S v we find. 



This can be also be expressed in terms of our standard 
null normals of Eq. as 



(54) 



-V^kfj, and C = VH^. 



with coefficients B 

Now we determine the rescaled null normals and fc M 
appearing in Eq. (|4"5)l . Since must be the same vector 
whether it is written in terms of and or and fc M , 
we have the relations 



B = B/f and C = fC, 

which together with Eq. (|46[) gives 

f _B_C__ 2 rfM H 
■' ~ B~ C ~ C dv ■ 



(55) 



(56) 



Evaluating the scale factor / requires knowledge of 
dMn/dv. It is straightforward to show that the area 
element ifq of S v changes along H as 

CvVn = -Ce {k)y /q, (57) 

so the change in the cross-sectional area An along H is 

dA H 



dv 



= - f ce [k) v^d 2 



(58) 



From the definition Mh = -J j4h/167t, it then follows that 
dM H 1 cL4 H 



dv 



s/MttA^ dv 



(59) 



The terms in the dynamical horizon flux law (|47[) arc 
calculated by noting that under the rescaling of the null 
normals and k^, 



and 



rf M In /. 



(60) 



The results are shown in Fig. [T^] from t ~ 10Mj to 
t = 20Mi. The energy flux of matter is neglected since 



we have T,, 



0. 



The flux associated with Ba^a 



labeled as "er^" flux" , is always the larger contribution to 
the growth of Mh, which is expected from the interpreta- 
tion of this term as a flux of gravitational energy. This is 
most pronounced along MTT2 and MTT4, with decreas- 
ing Mh during the evolution, and clearly indicates that 
their appearance is a consequence of the sufficiently high 
gravitational energy flux across them. We have seen in 
Sec. IV Bl that for weak gravitational waves and with the 
same gauge condition for the evolution, no such MTTs 
appear. The maximum number of MTTs that can exist 
at the same time may also be linked to the structure of 
the gravitational waves, as shown in the inset of Fig. [51 
although we have not explored this aspect further. 

The fluxes increase rapidly near each bifurcation point. 
This is because of our choice of normalization for in 
Eq. (HI]), which propagates into V^. To understand this, 
let us write as the spacctimc coordinates of S c that 
bifurcates, with foliation parameter v = c say. Then on 
a nearby S v , we can approximate dx l v / dt by 



dxl 
dt 



0_ 

dt 



A 

2 VTt 



1 



t c \ 



(61) 



where A is some function. As t — > t c , this quantity di- 
verges as does the norm of V 1 *, and leads to the higher 
values of the fluxes measured along V 11 . This singular 
behavior could be absorbed into a redefined foliation pa- 
rameter v' = v'(v). Also, any visible discontinuities in 
the fluxes across different sections of H in Fig. [12] are 
due to the difficulty in finding S c exactly (as indicated 
by the data points in Fig. [SJ even searching for MOTSs 
at every At = 0.01 is insufficient for this purpose). 



E. Angular Momentum Flux Law 



The angular momentum J defined in Eq. (|48|) depends 
on a choice of rotation vector </> M on S v . If S v is ax- 
isymmetric, the natural choice of is the axial Killing 
vector. In general spacetimes no such Killing vector ex- 
ists, but one can nevertheless define a suitable ^ [62j by 
requiring it to have closed orbits, and be divergence-free 



0. 



(62) 



This notion has been further refined to calculate approx- 
imate Killing vectors [2!| |3C| in black hole simulations, 
and we will make use of this choice here. They were also 
used to compute J of the initial data sets in Sec. [TTJ 
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Gourgoulhon has generalized the Damour-Navier- 
Stokes equation for null hypersurfaces to trapping hori- 
zons and used it to derive a flux law for the change in 
J along a hypersurface T~L foliated by 2-surfaces S v (not 
necessarily MOTSs) with foliation parameter v [2l| , 



dJ_ 

dv 



= - [ T, LV ^T v ^d 2 x 
Js v 

167r Js v 

= - / T^t" Jdd 2 x 
Js v 



(63) 



(64) 



1 
1 



S v 



->7T 



y/qd x 

[e^d^c-^Cv^] V^d 2 x, 



where the vectors V 1 = Bl* - Cfc' 1 and = Bl^ + 
Ck^ arc tangent and normal to %, respectively. The 
first integral in Eq. (|6"4"]) is the angular momentum flux 
due to matter. The second integral can be thought of as 
the flux due to gravitational radiation and vanishes if S v 
is axisymmetric. In addition, it is usually required that 
(j)^ be Lie transported along the dynamical horizon, 



(65) 



so that the last integral in Eq. (|64|) vanishes when S v is 
an MOTS [2l| . This requirement ensures that in the ab- 
sence of matter and gravitational radiation, the angular 
momentum flux will be zero along an MTT as expected, 
instead of there being some physically unmeaningful flux 
simply due to measuring J about different axes. 

Here we evaluate the angular momentum flux law for 
the dynamical horizon H found in Sec. IV CI for A = 0.5. 
Because we calculate J with </> M being an approximate 
Killing vector, Eq. (|65j) is not satisfied in general, and so 
we must keep the last integral in Eq. (|64|) . We use the 
same tangent vector and foliation parameter v along 
each section of H. as in Sec. IV Dl and the null normals to 
S v given in Eq. ([29]) . The values of the terms in Eq. (f64|) 
arc shown in Fig. [H from t = 10M; to t = 20M ; . The 
first integral is neglected since T^ v = 0. The two terms 

in the second integral are labeled as "Ba^a^ flux" and 
u Ca (k) a (<p) flux „ Thc lagt integral is labeled as 

flux" . The angular momentum flux dJ j dv is dominated 
by thc flux associated with Ba^la^^ , due to the large 
cr^l produced by the gravitational waves. Thc magni- 
tude of dJ/dv vanishes initially, becomes largest along 
the end of MTT1 and the beginning of MTT2 when the 
gravitational waves reach the black hole, and settles back 
down to zero again along the successive MTT sections. 
Because dJ/dv alternates sign along H, the net change 
in J turns out to be small. The terms in the angular 
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FIG. 13. Terms in the angular momentum flux law of Eq. (|64[) 
plotted against the foliation parameter v along each section 
of H. Along MTT1, we choose v = t. Along the other MTT 
sections, we choose v — on the first S v we find. 



momentum flux law also diverge near each S v that bifur- 
cates into two MTTs, just like the terms in the dynamical 
horizon flux law in Fig. 1121 and again is a consequence of 
our choice of as discussed at the end of Sec. IV Dl 



VI. THE EVENT HORIZON 



A. Basic Definitions and Concepts 



The standard definition of thc surface of a black hole 
is the event horizon, the boundary of the set of all points 
that are not in the causal past of future null infinity ■ 
It is a null hypersurface, generated by null geodesies that 
have no future endpoints. As defined, the event horizon 
is a 3-surface, but it is common to refer to the intersec- 
tion of this surface with St as the event horizon as well. 
In contrast to an MOTS, the event horizon can only be 
found after the entire future history of the spacetime is 
known. Because of its teleological nature, the event hori- 
zon can behave nonintuitively. For instance, before a 
gravitational collapse has occurred an event horizon al- 
ready forms, even though there is no flux of energy or an- 
gular momentum across it yet. In this section we describe 
our method of finding the event horizon, and contrast its 
properties with those of the MTTs found in Sec. [V] 



13 




14 14.4 14.8 15.2 

1_J 



25 50 75 100 125 

tIM. 

i 

FIG. 14. Irreducible masses of the event horizon Meh and 
the MOTSs Mh during the evolution with A = 0.5. At the 
very beginning of the evolution Meh is already increasing, 
while Mh is still fairly constant. As the inset shows, Meh 
grows very slightly when Mh changes the most. 

B. Event Horizon Finder 

The event horizon is located in a spacetime by follow- 
ing geodesies backward in time. It is well known |63| 
that null outgoing geodesies in the vicinity of the event 
horizon, when followed backwards in time, will converge 
onto the event horizon exponentially. Therefore, given a 
well-chosen congruence of geodesies, one can trace the 
event horizon of the spacetime with exponentially (in 
time) improving accuracy. 

Our event horizon finder [64[ tracks a set of geodesies 
backwards in time. The initial guess for the event hori- 
zon is chosen at some late time when the black hole is in 
a quasistationary state. At this time, the apparent hori- 
zon and event horizon coincide closely, and the apparent 
horizon is used as the initial guess. The initial direction 
of the geodesies is chosen to be normal to the apparent 
horizon surface, and the geodesies are integrated back- 
wards in time. The geodesic equation requires values for 
the metric and its derivatives for each geodesic at each 
point in time. These values are obtained by interpolation 
from the values computed during the evolution. With an 
appropriate form of the geodesic equation, we can follow 
a geodesic as a function of coordinate time t, rather than 
the affinc parameter along the geodesic. 



tude A — 0.5. The surface area Aeh of the event horizon 
is computed by integrating the metric induced on its sur- 
face by the spatial metric ■ The irreducible mass of the 
event horizon is then given as Meh = n/^eh/IBtt. This 
is shown in Fig. 1141 together with the irreducible mass 
Mh along the MTTs. An obvious difference is that Meh 
always increases in time, and the event horizon does not 
bifurcate like the MTTs shortly after t = UM t . The 
event horizon is also already growing at the very begin- 
ning of the evolution, before the gravitational waves have 
hit the black hole. By t = 14Mj, the value of Meh has 
almost doubled while Mh is still fairly close to its initial 
value. In fact, during the time when multiple MTTs are 
present and one would intuitively expect the black hole 
to be the most distorted, the event horizon shows very 
little growth. 

This peculiar behavior of the event horizon was also 
illustrated in [l8| for the gravitational collapse of spheri- 
cal dust shells, and explained with the null Raychaudhuri 
equation (l7j . 

* = -^ ) -^) CT (^-8 7 r V r j (66) 

where A is an affinc parameter along the congruence of 
null geodesies that generate the event horizon, with tan- 
gent vector I 11 . The area element y/h of the event horizon 
is related to the expansion 9{i) by dy/h/dX = On^s/h, and 
substituting this into Eq. (|66p gives 

^ = (^) - S-V) Vh. (67) 

In dynamical situations we will generally have 9^ =/= on 
the event horizon, and this accounts for its accelerated 
growth, which is evident even at early times in our evolu- 
tion when the shear a^l is negligible. When the pulse of 
gravitational waves hits the black hole, a^l on the event 
horizon becomes large, and according to Eq. (|67|) this will 
decelerate its growth, even causing the growth to become 
very small in our case. 

At late times, the event and apparent horizons eventu- 
ally coincide as both a^l and 9^ go to zero on the event 
horizon while the apparent horizon becomes null. Finally, 
Fig. [TS] shows a spacetime diagram of the event horizon 
and the dynamical horizon T-L, with the spatial dimension 
along the z— direction suppressed. The null generators of 
the event horizon are shown as dotted red lines, and lie 
outside the solid grey surface of T~L , except when they co- 
incide at late times. In Fig. [15] the event horizon's cross 
section appears to be shrinking at late times. The con- 
stancy of the area of the event horizon (cf. Fig. fTl)) shows 
that this is merely a coordinate effect. 



C. Contrasting the Event Horizon with MTTs 



VII. DISCUSSION 



We find the event horizon for the evolution in which 
the ingoing gravitational waves have the largest ampli- 



In this paper, we investigate marginally trapped tubes 
and the event horizon for rotating black holes distorted 
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FIG. 15. Spacetime diagram of the event horizon and dynam- 
ical horizons for A = 0.5. The dotted red lines are the null 
generators of the event horizon, while the solid grey surface 
represents the dynamical horizons. 



by a pulse of ingoing gravitational waves. For small dis- 
tortions (low amplitude A), the simulations do not ex- 
hibit any unexpected behavior: the area of the apparent 
horizon is initially approximately constant, it grows when 
the gravitational radiation reaches the black hole, and 
then settles down to a constant value after the highly 
dynamical regime is over. However, for strong distor- 
tions, we find much more interesting behaviors of the 
MOTSs. A new pair of MOTSs appears outside the origi- 
nal MOTS. These new surfaces are initially close together 
and move rapidly away from each other, indicating that 
at the critical time when they first appear they are coin- 
cident (although this particular event cannot be resolved 
in an evolution with finite time step). The inner surface 
of such a pair shrinks, eventually approaches the origi- 
nal MOTS, and then these two surfaces annihilate each 
other. For amplitude A = 0.4 this process happens once, 
for A = 0.5 this happens twice, and there is a short time 
interval during which five MOTSs are present in the sim- 
ulation. 

The MTTs traced out by the MOTSs are smooth, 
and appear to combine into one smooth hypersurface 
(although the critical points where different marginally 
trapped tubes combine with each other cannot be re- 
solved). When the black hole is distorted, we find that 
this hypersurface is everywhere spacelikc and a dynami- 
cal horizon. We investigate how the black hole grows by 
evaluating the dynamical horizon flux law of Ashtekar 
and Krishnan [IE], and find that the gravitational 
energy flux is largest across the sections of the dynami- 
cal horizon that decrease in cross-sectional area with in- 
creasing time. We also evaluate the angular momentum 
flux law of Gourgoulhon [2l[ along the dynamical hori- 



zon, but instead of using a rotation vector ^ that is 
Lie transported along the dynamical horizon, we use an 
approximate Killing vector [29[, since we prefer to cal- 
culate the angular momentum itself in this way. The 
angular momentum flux law is based on the general- 
ized Damour-Navier-Stokcs equation, which treats the 
black hole as a viscous fluid. Evaluating the generalized 
Damour-Navier-Stokes equation itself could aid in devel- 
oping physical intuition about black holes in numerical 
spacetimes. 

In illustrating the procedure for finding multiple 
MOTSs, caution must be taken to locate the apparent 
horizon with MOTS finders when the MOTS found at a 
previous time is used as an initial guess. If the MOTS 
finder is not run frequently enough, new MOTSs will be 
missed and an erroneous apparent horizon will be iden- 
tified. This raises the issue of whether the true appar- 
ent horizon was indeed located in similar work involving 
highly distorted black holes in the past (e.g. Q). A better 
understanding of the slicing dependence of the MOTSs 
in our simulations would also be helpful in choosing a 
more natural slicing condition that gives a single dynam- 
ical horizon that only grows in the cross-sectional area 
with time in highly dynamical situations. 

When computing the event horizon, we find it to be 
smooth, and enveloping the complicated structure of the 
MOTSs. As can be seen in Figs. [14] and [15] the event 
horizon is very close to the apparent horizon at late times, 
as one would expect. The motion of the event horizon is 
restricted by the fact that it is foliated by null geodesies. 
Therefore, in order to encompass the MOTSs, the event 
horizon begins to grow much earlier, and even at the start 
of our simulation the event horizon is already consider- 
ably larger than the apparent horizon. At early times, 
t < 10Mj, the event horizon approaches the apparent 
horizon exponentially. The rate of approach should be 
given by the surface gravity of the initial black hole, but 
we have not verified this in detail, as our simulation docs 
not reach sufficiently far into the past. This could be 
checked by placing the initial pulse of gravitational radi- 
ation at a larger distance from the black hole. The growth 
of the event horizon is described by the Hawking-Hartlc 
formula [65| , which may also be evaluated to give a more 
complete comparison of MTTs and the event horizon. 

Our findings are analogous to the behavior of MOTSs 
and event horizons in the Vaidya spacetime, as worked 
out in detail in the Appendix. In particular, for strong 
accretion, the Vaidya spacetime can also exhibit multiple 
MOTSs at the same time, all of which foliate dynamical 
horizons. Both in the Vaidya spacetime and our distorted 
Kerr spacetimes, the event horizon begins to grow much 
earlier before multiple MOTSs appear. By choosing a 
mass functions m(v) that has two strong pulses of ac- 
cretion, the Vaidya example in the Appendix would also 
produce five concentric MOTSs similar to that seen in 

Fig. m 
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Appendix: Multiple Horizons in the Vaidya 
Spacetime 

The ingoing Vaidya spacetime is a spherically sym- 
metric spacetime describing a black hole that accretes 
null dust [£f|. It shares similar features to the distorted 
Kerr spacetimes presented in this paper, which we men- 
tion here briefly. The ingoing Vaidya metric in ingoing 
Eddington-Finkelstein coordinates (v,r,8,(f>) is 



dg2 = ,_2mM\ ^ 



2dvdr 



? -dn 2 , (A.i) 



where v = t + r is advanced time (not to be confused 
with the foliation parameter v of dynamical horizon in 
the main text). From the Einstein equations, the stress- 
energy tensor is 



T — 



dm/dv 

-4^5-(M(^«)- 



(A.2) 



With the choice of radial outgoing and ingoing null vec- 
tors 



V 



1/ 2^ 
' 2 V r 



,0,0 



and = (0,-1,0,0) 



normalized so that l^k^ 
normals arc 



(A.3) 

T, the expansions of the null 



'(0 



= - 1 



2m(v) 



and 6t k ) = — . 

r 



(A.4) 



From this, we see that MOTSs are located at r 
2m(v), or 



(A.5) 



The number of solutions to Eq. (|A.5[) . i.e. the number of 
MOTSs, can be conveniently discussed with the diagram 
shown in Fig. 1161 The thick solid lines represent three 




FIG. 16. Mass functions m(v) of the Vaidya spacetime for 
three amplitudes A — 0.25, 0.5, and 1, along with the straight 
lines (v — t)/2. MOTSs exist at the intersections of these func- 
tions. For A = 0.5 and 1, there are up to three intersections, 
as illustrated by the dashed black line which intersects the 
A = 1 mass curve three times. 



different mass functions m(v) plotted vs v. The right- 
hand side of Eq. (|A.5[) is a family of straight lines (one for 
each t) represented by the thin diagonal lines in Fig. 1161 
For a given t, the number of intersections between the 
(v — t)/2 and the m(v) curve gives the number of MOTSs 
at that particular t. The straight line |(w — i) has slope 
1/2, so if dm/dv < 1/2 for all v, then there will be exactly 
one intersection^ for every t. If 



dm 1 
~dv > 2 



for some v, 



(A.6) 



then the m(v) curve will have regions that are steeper 
than the straight line. By adjusting the vertical inter- 
cept of the straight line, equivalent to choosing a suit- 
able t, the straight line will pass through a point with 
dm/dv > 1/2. At this point, m(v) passes from below to 
above the straight line, so there must be an additional 
intersection at both smaller and larger v, for a total of 
three MOTSs. Thus, sufficiently rapid mass accretion 
(large dm/dv) results in multiple MOTSs. 

The signature of a spherically symmetric MTT de- 
pends on the sign of [l!| 



C 



rp mil/ 



l/^A^-T^lW 



(A.7) 



1 Assuming m(v) is non-decreasing, and has finite bounds for v 
±oo. 
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FIG. 17. Locations of MOTSs (solid lines) and event horizons 
(dashed lines) in the Vaidya spacetime. For A = 0.25 there 
is one MOTS at all times. For A — 0.5 and 1, up to three 
MOTSs exist at a time t. The event horizons approach the 
MTTs at very early and late times, and start growing much 
earlier than the MTTs. The inset shows a larger interval in t. 



where An is the cross-sectional area of the MTT. The 
MTT is spacelike if C > 0, null if C = 0, and timelike if 
C < 0. From Eq. (|A3j) and Eq. (fO) . 



T lu/ m v = and T^Pk" = 0, (A. 



so we see that C > for the Vaidya spacetime as long 
as dm/dv > 0. Furthermore, since 9^) < 0, these MTTs 
will also be dynamical horizons. 

The event horizon is generated by radial outgoing null 



geodesies satisfying 

dr _ 1 
dv ~ 2 



2m(v) 



(A.9) 



Integrating this differential equation requires knowledge 
of the event horizon location at some point. This is usu- 
ally supplied by the final state of the black hole, when 
accretion has ended. 

To close, we illustrate these considerations with a con- 
crete example. We choose the mass function 



m(v) 



m , 



mo 



Atrial) 



W 2 ' 



v < 



v > 



(A.10) 



similar to that presented in [67[ {Am is the mass ac- 
creted by the black hole, and W determines the time scale 
of accretion). We set m = 1, W = 0.5, and consider 
three different amplitudes A = 0.25, 0.5, and 1. FigurefTBI 
shows the respective mass functions, and we see that 
A = 0.25 never leads to multiple MOTSs, while A = 1 
clearly exhibits three MOTSs for certain t. It is easy 
to show that Eq. (TA~6]) implies Am > 4W A /(3v / 3). The 
locations of the MOTSs in (r, t) coordinates are shown 
in Fig. OH For A = 0.25, there is only one MOTS at 
all times. For A = 0.5, there are up to three MOTSs 
at a single time. A new MOTS appears at r = 2.5 im- 
mediately after t = — 2, and bifurcates into two MTTs. 
One of these MTTs shrinks and annihilates with the in- 
nermost MTT at t = —1.93256, while only the outermost 
MTT remains at late times and grows towards r = 3. For 
A = 1, there are again up to three MOTSs at a single 
time, but a new MOTS appears earlier at t = —2.63822. 
After t = —1.96824, only one MOTS remains and grows 
towards r = 4. Also shown in Fig. [T7] are lines of con- 
stant v indicating when accretion begins (v = 0), and 
when m(v) has increased by 50% and 80%, respectively 
(v = W and v = 2W). 

The event horizons for the three cases are computed 
by integrating Eq. (|A.9|) backward in time, starting with 
rEn(v — > oo) = 2(1 + A) mo. The resulting surfaces are 
shown as the dashed curves in Fig. [T7J The event horizon 
is located at r = 2 in the far past, starts growing long 
before m(v) increases, and asymptotically approaches the 
MTT of the final black hole for all amplitudes A. 
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